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Abstract 

Correlated ab initio electronic structure calculations are reported for the 
polymers lithium hydride chain [LiH]^ and beryllium hydride [-E^-f^jco- 
First, employing a Wannier-function-based approach, the systems are stud- 
ied at the Hartree-Fock level, by considering chains, simulating the infinite 
polymers. Subsequently, for the model system \LiH\oa, the correlation effects 
are computed by considering virtual excitations from the occupied Hartree- 
Fock Wannier functions of the infinite chain into the complementary space 
of localized unoccupied orbitals, employing a full-configuration-interaction 
scheme. For [-Be2-ff4]ocM however, the electron correlation contributions to 
its ground state energy are calculated by considering finite clusters of increas- 
ing size modelling the system. Methods such as M0ller-Plesset second-order 
perturbation theory and coupled-cluster singles, doubles and triples level of 
theory were employed. Equilibrium geometry, cohesive energy and polymer- 
ization energy are presented for both polymers, and the rapid convergence 
of electron correlation effects, when based upon a localized orbital scheme, is 
demonstrated. 

I. INTRODUCTION 



Polymers represent a class of one-dimensional infinite crystalline systems where ab initio 
Hartree-Fock (HF) self-consistent field (SCF) methods are well developed^. An available 
program package is CRYSTAL!. However, in order to be able to calculate the structural 
and electronic properties of polymers with an accuracy that allows a meaningful comparison 
with experiment, it is usually necessary to include the effects of electron correlations into the 
theory. The most widely used approach here is density-functional theory (DFT). Despite 
its indisputable success in solid state physics and computational chemistry as a computa- 
tionally cheap routine tool for large-scale investigations, DFT has the drawback that results 
depend highly on the chosen functional, and cannot be improved in a systematic way. Wave- 
function-based quantum chemical ab initio techniques on the other hand are free from this 
flaw, and provide a large array of methods of different accuracy and computational cost. 
Thus it is desirable to extend their applicability to infinite systems such as polymers. 
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Electron correlations are mostly a local effect and therefore localized molecular orbitals 
are preferable to the canonical HF solutions for the treatment of large molecules.! Similarly, 
in infinite systems (localized) Wannier functions provide a better starting point for an ab 
initio treatment of electronic correlations than the (canonical) Bloch functions. Previous 
studies of polymers obtained the Wannier orbitals from an a posteriori localization of the 
Bloch functions according to a given prescription.! During the last years, in our group a HF 
approach was developed which allows the direct determination of Wannier orbitals within 
the SCF process.! Various applications to one- and three-dimensional infinite systems proved 
the numerical equivalence of our Wannier-function-based HF approach to the conventional 
B loch-function-based counterpart .!'! 

In this paper HF-SCF calculations and subsequent correlation energy calculations 
are presented for the lithium hydride chain [LiH]^ and the beryllium hydride polymer 
[Be 2 H i ] 00 . As a simple, but due to its ionic character, nontrivial model polymer, the 
lithium hydride chain system has been previously dealt with in a number of studies.i~!. 
In the present contribution we extend our previous calculation ! to a wave-function-based 
ab initio study of electron correlation effects using a combination of the full configuration 
interaction (FCI) method and the the so-called incremental scheme.11^0 The latter approach 
consists basically in an expansion of the total correlation energy per unit cell in terms 
of interactions of increasing complexity among the electrons assigned to localized orbitals 
(Wannier functions) comprising the polymer under consideration. The electron correlation 
energy increments needed to establish the total energy per unit cell are evaluated by con- 
sidering virtual excitations from a small region of space in and around the reference cell, 
keeping the electrons of the rest of the crystal frozen at the Hartree-Fock (HF) level. The 
fast convergence of the incremental expansion allows to truncate it at relatively low order 
and thus to calculate the correlation energy of an infinite system without modelling it as 
a finite cluster. However, neither the FCI method nor the incremental approach based on 
polymer Wannier orbitals can at present be used for systems with a more complicated unit 
cell. Therefore, the second system investigated by us, the beryllium hydride polymer, was 
treated at the coupled-cluster (CC) and M0ller-Plesset second-order perturbation (MP2) 
level of theory. Starting from the Wannier HF data the correlation corrections to the total 
energy per unit cell were derived from quantum chemical calculations of finite model systems 
using the MOLPRO molecular orbital ab initio program package.0 To our knowledge, this 
system was studied at the HF level two decades ago by Karpfen using the crystal orbital 
method, i.e., without including correlation effects.t3 Recently, its monomer beryllium dihy- 
dride BeH 2 has been well characterized theoretically using reliable ab initio and density 
functional theory methods-Hl! 

The remainder of the paper is organized as follows. In section |H| the applied methods 
are briefly described. The calculations and results are then presented in section |HI[ Finally, 
a summary is given in section [TV]. 



II. APPLIED METHODS 



Section [11 A| gives brief outline for the theory within a restricted HF (RHF) framework. 
Sections [II B| and |II Cj , respectively, describe the incremental scheme and a simple approach, 
to compute electron correlation effects in polymers. 
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A. Wannier orbital based Hartree Fock approach 



Our approach, described in more detail in previous publications.oH is based upon the 
direct determination of the orthonormal Wannier-type (localized) orbitals for the polymer. 
Denoting by | a(Rj)) the Wannier orbitals of a unit cell located at lattice vector Rj, the 
set {|a(Rj));a = l,n c ;j = l,N} spans the occupied HF space. Here, n c is the number of 
orbitals per unit cell, and iV(— > oo) is the total number of unit cells in the system. In our 
previous work we showed that one can obtain n c RHF Wannier functions, {\a), a = l,n c } 
occupied bv 2n c electrons localized in the reference unit cell (denoted C) by solving the 
equations@~! 



(T + U + YpJp -^)+EE A^| 7 (R fe ))(7(R fe )|)|a) = e», W 

13 kGAf 7 

where T represents the kinetic-energy operator, U represents the interaction of the electrons 
of C with the nuclei of the whole of the crystal, while Jp, Kp, respectively, represent the 
Coulomb and exchange interactions felt by the electrons occupying the (3-th Wannier function 
of C, due to the rest of the electrons of the infinite system. The first three terms of Eq.(|l]) 
constitute the canonical Hartree-Fock operator, while the last term is a projection operator 
which makes the orbitals localized in C orthogonal to those localized in the unit cells in 
the immediate neighborhood of C by means of infinitely high shift parameters A^'s. These 
neighborhood unit cells, whose origins are labeled by lattice vectors R^, are collectively 
referred to as N ' . The projection operators along with the shift parameters play the role 
of a localizing potential in the Fock matrix, and once self-consistency has been achieved, 
the occupied eigenvectors of Eq.([I]) are localized in C, and are orthogonal to the orbitals 
of M — thus making them Wannier functionsEHl As far as the orthogonality of the orbitals 
of C to those contained in unit cells beyond Af is concerned, it should be automatic for 
systems with a band gap once M has been chosen to be large enough. Based upon our past 
experience regarding a suitable choice of A/",^i in the present calculation we included up 
to the third nearest-neighbor unit cells in N '. For the details concerning the computation 
of various terms involving lattice sums (U, J, and K) involved in Eq. (P for the case of 
polymers, we refer the reader to reference!. 



B. Incremental method 

Electron correlation effects in the ground states of a large number of three-dimensional 
ionic and covalent solids,@ as well as polymers^ have been studied with the incremental 
scheme. All these calculations used localized orbitals of finite clusters as a basis set for the 
correlation treatment. In the present work on the lithium hydride chain we use directly the 
Wannier-functions of the infinite system. A related study of the three-dimensional lithium 
hydride solid has been published elsewhere0. 

The correlation energy per unit cell is expanded as 

E corr = ^ei+ ^2 Ae lj + Ae ijk + ... (2) 

i <ij> <ijk> 



3 



where the summation over i involves Wannier functions located the reference cell, while those 
over j and k include all the Wannier functions of the crystal. The "one-body" increments £j 
= A^j are computed by considering virtual excitations only from the z-th Wannier function, 
freezing the rest of the polymer at the HF level. The "two-body" increments Ae^ are defined 
as Aeij = Eij — (Affj + Asj) where is the correlation energy of the system obtained by 
correlating two distinct Wannier functions i and j. Thus Ae^ represents the correlation 
contribution of electrons localized on two "bodies" i and j. Higher-order increments are 
defined in an analogous way. Finally, summing up all increments, with the proper weight 
factors (according to their occurrence in the unit cell of the polymer), one obtains the 
exact correlation energy per unit cell of the infinite system. In order to get reliable results 
a size-extensive correlation method should be used, although non size-extensive schemes 
also may provide reasonable estimates if the incremental expansion is truncated at low 
order. In the present work for the lithium hydride chain we choose the strictly size-extensive 
full configuration interaction (FCI) method. As mentioned earlier, when computing the 
correlation contributions via Eq. (|2|), except for the orbitals involved (say orbitals % and j 
for the two-body increment Ae^), the rest of the occupied Wannier orbitals of the infinite 
solid are held frozen at the HF level. The region containing these frozen orbitals plays the 
role of the "environment" for the electrons involved in the correlated calculations, and its 
contribution can be absorbed in the so-called "environment potential" JJ env defined as 

U£"= E (2(pa(R,)|-^|ga(R J ))-(p«(R,)|-^|a(R J )g)), (3) 

where £ represents the unit cells of the environment, p and q are two arbitrary basis functions, 
and the factor of two in the first term is due to the spin summation. The sum of Eq.(|3]) 
involves infinite lattice sum over the environment unit cells, and is computed by simply 
subtracting from the lattice summed J and K integrals (cf. Eq. ([j])) obtained at the end 
of the HF iterations, the contributions corresponding to the orbitals being correlated. Once 
Up^ Y has been computed, one is left with an effective Hamiltonian involving a finite number 
of electrons located in the region whose Wannier orbitals are being correlated. Physically 
speaking U^ Y represents the influence of the environment electrons on the electrons being 
correlated, explicitly. In the present calculations the Li Is 2 core shell was also kept frozen, 
and its contribution was also included in U^ v . The basis functions p and q were restricted 
to those of the reference cell and the adjacent cells up to the third- nearest neighbors. 

The virtual orbitals used for computing the correlation effects were also localized. They 
were obtained by first orthogonalizing the basis set to the occupied space by using corre- 
sponding projection operators, as suggested by PulayS Subsequently the basis functions 
are orthogonalized to each other using the symmetric-orthogonalization procedure, yielding 
a localized and orthonormal virtual orbital set .El The number of virtual orbitals per unit 
cell considered for a specific increment corresponds to the number of basis functions per 
unit cell minus the number of occupied orbitals per unit cell. The virtual orbitals have been 
expanded in the same basis set as described above for U^ Y . 
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C. A simple approach 



In principle the total energy E tot per [Be^H^ unit cell of beryllium hydride may be 
obtained as the limit 

E = lim E(jge2 " +lg4 " +2) , (4) 

i.e., by performing calculations for increasingly long oligomers H(BeH 2 )2 n BeH. In order to 
reduce finite-size effects due to the termination of the oligomers by one beryllium and two 
hydrogen atoms saturating the dangling bonds of ■{BeH2)2n , one may consider instead 



E = lim AE n = lim 

n— >oo n— »oo 



H 



4n+2) 



(5) 



i.e., the energy change between subsequent oligomers differing by a single unit cell. There- 
fore, identical unit cells were used as building blocks for both oligomers, i.e., the geometrical 
optimization was restricted only to parameters relevant for the polymer beryllium hydride. 
Since the convergence of AE n with respect to n is much faster for the correlation contri- 
butions than for the HF energy, and HF programs treating the infinite system are at hand 
(CRYSTAL, WANNIER), we use Eq. (H) only for the correlation energy per unit cell. This 
approach has previously been used successfully in calculations for trans-polyacetylene,0 and 
some boron-nitrogen polymers.ll! 



III. CALCULATIONS AND RESULTS 

A. [LiHU 

HF ground state calculations are a necessary prerequisite for the application of the in- 
cremental approach to electron correlation. We performed such calculations for a lithium 
hydride chain oriented along the x-axis using the WANNIER code0. The reference cell con- 
tained hydrogen at the (0,0,0) and lithium at the (a/2, 0, 0_L where a is lattice constant. We 
adopted the extended basis set optimized by Dovesi et alx3 First, all-electron Wannier HF 
calculations were performed at the different lattice constants in the range 2.8-4.0 (A) and 
the total HF energy per unit cell for various lattice constants near the equilibrium was fitted 
to a cubic polynomial in order to derive the ground state HF equilibrium lattice constant and 
total energy. After determining the Wannier orbitals for each value of the lattice constant, 
the corresponding FCI calculations were performed by means of the incremental scheme. 
The expansion of the correlation energy per unit cell was restricted to one- and two-body 
increments, and included interactions up to third-nearest neighbor unit cells. Contributions 
from higher order increments as well as from interactions between more distant cells proved 
to be negligible. The equilibrium values for the FCI energy per unit cell and the lattice 
constant were determined as described for the HF results. The main contribution of 98.8 % 
to the correlation energy per unit cell at the equilibrium geometry {E=— 0.0307a. u.) comes 
from the one-body term. Two-body terms for first-, second- and third-nearest neighbors 
contribute with 1.15, 0.01 and 0.001 %, respectively. Our results are summarized in table 
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0. It is quite obvious from table p] that, as a function of distance, the two-body correlation 
effects converge very rapidly. 

Since the Li basis set used here is suitable only for the ionic LiH molecule, we cannot 
get a good result for the atomic reference energy of the neutral Li atom (which is needed to 
determine the cohesive energy). Therefore, for this almost ideally ionic chain the cohesive 
energies both at the HF and the correlated level are obtained by subtracting the electron 
affinities (EA) and ionization potential (IP) from the dissociation energy calculated with 
respect to the ions Li + and H~ . The HF values of EA and IP are determined using the finite- 
difference atomic HF program MCHF§1. The experimental values of EA and IP were taken 
as the CI limit, i.e., disregarding the very small relativistic effects. For the polymerization 
energy we optimized the Li-H distance for the 1 S + ground state of the monomer at the 
HF and CI level. Our results are summarized in table JTI]. It is clear from table [TT] that, as 
expected, correlation effects contribute significantly to the cohesive energy. However, they 
do not make any significant contribution to the lattice constant of the system. 



B. [fierce 

Beryllium hydride has attracted considerable interest as a rocket fuel on account of its 
high heat of combustion. It has also been considered as a moderator for nuclear reactors. 
From the previous studiesii we also know that it is poisonous and difficult to prepare for 
experiment. Even though there is no or very little experimental information about the 
polymer.it has been studied theoretically using reliable ab initio methods at the HF level by 
Karpfenlii In the present work we have studied this polymer at the HF and the correlated 
level. The Wannier-orbital-based HF-SCF approach, coupled-cluster (CC), and M0ller- 
Plesset second-order perturbation (MP2) theory were employed to determine the equilibrium 
structures and total energies per unit cell. In our calculations the unit cell included two 
beryllium and four hydrogen atoms and has a perfect tetrahedral structure with all four Be- 
H bond distances equal, i.e., there are two HBeHBeH planes that are perpendicular to each 
other, with the beryllium atoms in their crossings. In the cluster approximation the unit 
cell is terminated by one beryllium and two hydrogen atoms. In this structure the terminal 
beryllium atoms have trigonal coordination while all others are distorted tetrahedrons. First 
we optimized the structure of this polymer at the HF-SCF level using the CRYSTALS 
program. The total HF energies obtained with the CRYSTAL program were then taken as 
an input for a re-optimization at the MP2, CCSD (CC singles and doubles) and CCSD(T) 
(CCSD with a perturbative estimate of triples) level. The correlation energy contributions at 
each geometry have been calculated with the MOLPRO molecular orbital ab initio program 
packaged by using the simplified finite-cluster approach in which we put n=3 in Eq. ([|). In 
this system the correlation energy converges rapidly with respect to cluster size, i.e., for n=3, 
one finds AE4 — A^^IO -6 a.u.. We have optimized the beryllium-hydride bond length 
( r Ben) and the lattice constant (a). We adopted polarized valence double-zeta ( 6-31G**) 
basis sets for beryllium and for hydrogen. The polarization functions consisted of a single 
p-type exponent of 0.75 Bohr -2 on hydrogen and single d-type exponents of 0.4 Bohr -2 
on beryllium. In our HF calculations for polymers we optimized the most diffuse s-type 
exponent, which is less than 0.1 in the original 6-31G** basis set, and obtained 0.15. A 
smaller value causes linear dependencies in the basis set when applied in the infinite system. 
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We have also calculated the cohesive energy per unit cell at the HF and correlated level. 
The atomic HF-SCF, MP2, CCSD and CCSD(T) reference energies (Be: -14.5668 a.u., 
-14.5928 a.u., -14.6131 a.u. and -14.6131 a.u.; H: -0.4982 a.u.) were obtained with the 
original 6-31G** basis sets. In addition to the cohesive energy, we have also calculated the 
polymerization energy. The geometry of the monomers was optimized at the SCF, MP2, 
CCSD, and CCSD(T) levels of theory employing the MOLPRO programEl. Our final results 
are summarized in table [TTTj . Due to the absence of experimental data or theoretical results 
at the correlated level, we compare our result only at the HF level. To the best of our 
knowledge, only Karpfenll! has performed a geometry optimization for this polymer within 
an ab initio crystal Hartree-Fock approach and his results are also given in table [TTT| . Our 
beryllium-hydrogen bond length is in good agreement with the one obtained by Karpfen, 
but our HF energy is lower 0.05 a.u. than the value of KarpfenEl. A possible reason is the 
use of d functions in our basis sets. 

IV. SUMMARY 

In conclusion, given a well-localized basis set of Wannier orbitals size-extensive standard 
quantum chemical methods such as full configuration interaction, coupled-cluster or many- 
body perturbation theory can be applied to evaluate ground state properties of polymers. 
Rapid convergence of the incremental expansion of the correlation energy is obtained for 
ionic systems, e.g., the simple model of the lithium hydride chain. In beryllium hydride 
polymer electron correlation accounts for 12-14% of the cohesive energy and 22-24% of 
the polymerization energy at all three levels of theory and reduces the lattice constant. 
In all the cases it was demonstrated that the use of localized orbitals leads to a rapid 
convergence of electron correlation effects, thus making it possible for one to compute the 
electron correlation effects of infinite systems. 
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TABLES 



TABLE I. Various increments to the correlation energy (in Hartrees) computed by the Wan- 
nier-function-based approach presented in this work. The results refer to lattice constant of 3.30 
(A). NN stands for nearest neighbors. 



Correlation 


Energy 


Increment 




one-body 


-0.0303345 


two-body (INN) 


-0.0003538 


two-body (2NN) 


-0.0000035 


two-body (3NN) 


-0.0000003 



TABLE II. Total energy Etot (Hartree), cohesive energy AE^ (eV), polymerization energy 
AE po ; (eV) per unit cell and lattice constant a (A) of the lithium hydride chain. 

Method E tot AE coh AE pol a 

WANNIER SCF -8.038047 3.8760 1.8067 3.3273 

CRYSTAL SCF -8.038031 3.8759 1.8063 3.3274 

FCI -8.068744 4.6545 1.4854 3.3300 



TABLE III. Total energy Etot (Hartree), cohesive energy AE co ^ (eV), polymerization energy 
AEp i (eV) per unit Be^H^ and lattice constant a (A), Be-H distance h (A) of beryllium hydride. 



Method 


Etot 


AE co /j 


AE po ; 


a 


h 


CRYSTAL SCF 


-31.6300 


13.70 


2.645 


3.958 


1.467 


MP2 a 


-31.7608 


15.85 


3.445 


3.958 


1.456 


CCSD a 


-31.7908 


15.56 


3.402 


3.969 


1.457 


CCSD(T) a 


-31.7944 


15.66 


3.478 


3.968 


1.458 


Karpfen 6 


-31.5780 






4.024 


1.470 



a correlation contributions added to CRYSTAL SCF energies. 

b performed with 7, 1/4 basis sets considering third neighbor's interactions. 
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